NUMERICAL SOLUTION OF THE SMALL DISPERSION LIMIT OF 
THE CAMASSA-HOLM AND WHITHAM EQUATIONS AND 
MULTISCALE EXPANSIONS 



S. ABENDA*, T.GRAVAt, AND C. KLEIN* 

Abstract. The small dispersion limit of solutions to the Camassa-Holm (CH) equation is 
characterized by the appearance of a zone of rapid modulated oscillations. An asymptotic description 
of these oscillations is given, for short times, by the one-phase solution to the CH equation, where the 
branch points of the corresponding elliptic curve depend on the physical coordinates via the Whitham 
equations. We present a conjecture for the phase of the asymptotic solution. A numerical study of 
this limit for smooth hump-like initial data provides strong evidence for the validity of this conjecture. 
We present a quantitative numerical comparison between the CH and the asymptotic solution. The 
dependence on the small dispersion parameter e is studied in the interior and at the boundaries of 
the Whitham zone. In the interior of the zone, the difference between CH and asymptotic solution is 
of the order e, at the trailing edge of the order and at the leading edge of the order e^/^. For the 
latter we present a multiscale expansion which describes the amplitude of the oscillations in terms of 
the Hastings-McLeod solution of the Painleve II equation. We show numerically that this multiscale 
solution provides an enhanced asymptotic description near the leading edge. 
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1. Introduction. The Camassa-Holm (CH) equation 
ut + {3u + 2iy)ux - e^{uxxt + '2uxUxx + uuxxx) = 0, x e M, t > 0, (1.1) 

was discovered by Camassa and Holm [5l as a model for unidirectional propagation of 
waves in shallow water, u{x, t) representing the height of the free surface about a flat 
bottom, V being a constant related to the critical shallow water speed and e a constant 
proportional to the mean water depth [13| . Equation had been previously found 
by Fokas and Fuchssteiner [15j using the method of recursion operators and shown 
to be a bi-hamiltonian equation with an infinite number of conserved functionals. 
It was also rediscovered by Dai [9] as a model for nonlinear waves in cylindrical 
hyperelastic rods, with u{x,t) representing the radial stretch relative to a pre-stressed 
state. Equation (jl.ip finally also arises in the study of the motion of a non-Newtonian 
fluid of second grade in the limit when the viscosity tends to zero [Ij . A class of two- 
component generalizations of the CH equation has been recently obtained in [14j . 
The initial value problem for (|l.ip 

u{x, 0) ~ uq{x), a; e M, 

presents interesting features: first ioi v = there may exist peakons, i.e., non- 
smooth solutions, second even for a smooth initial datum uo{x) the wave-breaking 
phenomenon may occur, that is the solution u(x,t) remains bounded while its slope 
becomes unbounded in finite time. This phenomenon was first noticed for ly — by 
Camassa and Holm ^ who showed that for smooth and odd initial datum ^0(2^) such 
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that uq{x) > for x < and Uo(0) < 0, the slope Ux{x,t) is driven to — oo in finite 
time. In the case v ~ 0, under the hypothesis that 

moix) := uo{x) - e^dxxUo{x) 

is smooth and summable, McKean [23j proves that the wave-breaking phenomenon 
occurs if and only if some portion of the positive part of toq (x) lies to the left of some 
portion of the negative part of toq {x) ■ The wave breaking phenomenon occurs also for 
v ^ and, in particular Constantin [6 shows that for initial data uo{x) in the Sobolev 
space H^{R), s > 3 there exists a unique solution to u e C"([0, r[, i/"(R)) n 

C^{[0,T[,H'^~^{M)) defined for some maximal time T > 0; moreover T < oo if and 
only if liminff^T{min2:gR Ux{x, t)} — — oo, i.e., for v ^ singularities in the solution 
may arise only in the form of wave breaking. 

In this manuscript we are interested in studying the behaviour of the solution 
t, e) of the Cauchy problem of the CH equation as e — > for smooth initial data 
for which breaking does not occur. To this end, we suppose that Uf){x) is in the 
Schwartz class, has a single negative hump and satisfies the following non-breaking 
condition 

mo(a;) + i/>0, a; e K. (1.2) 

In this case for any <: > 0, m{x,t;e) — u{x,t) — e^dxxu{x^t) is in the Schwartz class 
and m(x, t]e) + v > Q [7j. 

For initial data in this class and e ^ 1, Grava and Klein [T7] show that the 
numerical solution u{x,t,e) of (|l.ip develops a zone of fast oscillations, as for the 
small dispersion limit of the Korteweg-de Vries (KdV) equation, see Fig. 11.11 

In this paper, following the works of Gurevich and Pitaevskii [5D], Lax and Lev- 
ermore |25| . Deift, Venakides and Zhou [10] . we claim and give numerical support for 
our claim that the description of the small dispersion limit of the CH equation follows 
l)-3) below: 

1) For Q < t < tc, where tc is a critical time, the solution u{x, t, e) of the CH 
Cauchy problem is approximated as e ^ 0, by u{x, t) which solves the Hopf equation 

+ (3u + 2i^)ua; = 0. (1.3) 

Here tc is the time when the first point of gradient catastrophe appears in the solution 
to the Hopf equation 

M(x,t) = uo(0, x = {-iuo{0 + 2u)t + i, 
and it is given by the relation 

1 

tc = . 

max5gK[-3u;)(0] 

2) There exists T > tc such that for tc < t < T , the solution of the CH equation 
is characterized by the appearance of an interval of rapid oscillations. If f'l'{u) < 
where /_ is the inverse of the decreasing part of the initial datum uo{x) and Uc is 
the critical point, these oscillations are described in the following way. As e ^ 0, the 
interval of the oscillatory zone is given by [x^ (t) , x^ (t)] where x^{t) are determined 
from the initial datum and satisfy the condition x^{tc) ~ x^{tc) = Xc, with Xc 
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Fig. 1.1. Solution to the CH equation for uq 
values of t. 
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the coordinate of the point of the gradient catastrophe of the Hopf solution. In the 
(x,t) plane x^if) describe a cusp-shape region. Outside the interval \x~{i), the 
leading order asymptotics is given by the solution of the Hopf equation. Inside the 
interval (x^ (t), (t)) we claim that the solution u{x, t, e) is approximately described, 
for small e by the one-phase solution of CH which may be expressed in implicit form 
in terms of elliptic functions as 



u{x, t, e) ~ f3i + (32 + [33 + - 2(/33 + i^) 
C = - (a; - (/3i /52 + /53 + 2iy)t - 



A(s,p) 
K{s) ^ 

q) — lixz 



1(3 dz ■ds,{z +p,t)'' 



fclog 



(1.4) 



with /3i > /32 > /Sa > —v, K(s) and A(s,p) the complete Jacobi elliptic integrals of 
the first and third kind of modulus s, respectively, where 



il3i + iy){l32 - p3) ^ P2-I33 

(/32 + i')(/3i-/33)' ^ (32 + iy' 



(1.5) 



The quantities Iq , p and the wave-number k are Abelian integrals defined respectively 
in ([231), ([23]) and (|2TT|) . the phase shift g(/ 3) is de fined in ((331). i?3(z, r) is the third 
Jacobi theta function of modulus t = iK{^/l — s'^)/K{s) (see [32] and p.8|) ). 

For constant values of the Pi , (II. 4p is an exact solution to CH in implicit form (see 
[2] and section [21). Moreover, for —v < < P2 < (3i, such solution is real-periodic 
and analytic in (, so it is appropriate to call it the one-phase solution. However, unlike 
the KdV case, it cannot be extended to meromorphic function on the C complex plane. 
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Indeed, the r.h.s of the first equation in (jl.4p is an eUiptic function in z, while C — (^(z) 
is real analytic and invertible for real z, but not meroniorphic in C. 

In the description of the leading order asymptotics of u{x,t,e), as e ^ 0, the 
quantities f3i depend on x,t and evolve according to the CH modulation equations 
which were derived in [T| for the one-phase periodic solution. In terms of the Riemann 
invariants Ps < (32 < Pi they take the form 

a*A + C,(/3)a,A = 0, C,;(/3) = ^^, z = l,2,3, (1.6) 

with/3= (/3i,/32,/?3) andcj = (^i +^2+/33+2i^)fc. Unlike the KdV case, the Whitham 
equations for CH are not strictly hyperbolic and this fact gives some technical diffi- 
culties in their numerical and analytical treatment. 

3) Near the left boundary of the cusp-shape region, in the double-scaling limit 
X x^ (t) and e — s- 0, in such a way that 

{x-x-{t))e-^^'^ 

remains finite, the asymptotic solution of the CH equation is given by 

u{x, t, e) = u{x^{t),t) + eia{x, t) co&{ip{x, t, e)/e) + 0{ei), 

where u{x, t) is the solution of the Hopf equation (|1.3p . the phase is given in (|4.5p 
and the function a{x, t, e) is, up to shifts and rescalings, the Hastings-Mcleod solution 
to the Painleve-II equation |2jy 

A,_^ = zA + 2A^, A,z€C, 

determined uniquely by the boundary conditions 



A{z) w Ai{z) {z +00 ), A{z) w y^-z/2 (z^-oo) 

with Ai{z) the Airy function. Such a solution is real and pole free for real values of z. 

We verify numerically the validity of the above asymptotic expansions 1), 2) and 
3) for the initial datum uo{x) = — 1/ cosh^ x and different values of v, see for instance 
Fig.O 

Remark. The asymptotic description of u{x, i, e) as e — ^ near the critical time 
t = tc has been conjectured in [TT] and studied numerically in [TH]. The asymptotic 
description of u{x,t,e) as e ^ at the right boundary of the cusp-shape region, 
namely near 2;+ {t) has not yet been studied even for the KdV case. 

The paper is organized as follows: In the next section we compute the one-phase 
solution to CH, in section [3] we obtain the small amplitude limit of the one-phase 
solution to the CH equation. In section [?] we perform a multiple scale analysis of 
the CH equations near the leading edge x~{t) of the oscillatory zone. In section [5] 
we present a numerical comparison of the small dispersion limit of the CH solution 
with the asymptotic formula (|1.4p . A quantitative numerical comparison of the CH 
solution and the multiscale solution is presented in section [6l 

2. The one-phase solution to the Camassa Holm equation. Let us look 
for a one-phase real-periodic travelling wave solution to (|l.ip of the form 

I ^\ t o \ o kx-ujt + (j)Q uj 

u{x,t) ^ {c-2iy) -2r](C), , c=- 

e k 



Fig. 1.2. Solution to the CH equation for the initial datum no = —sech}x and v = 1.2 at 
t = I in blue and the corresponding asymptotic solution via the Hopf equation and the one-phase 
CH solution solution. 



where k is the wave number, w the frequency and (po is a phase to be determined from 
the initial conditions. When we plug r]{Q into the CH equation (|l.ip . after integration 
we get 

k^{r] + i^)ril ^Tj^ + {c- 2v)rf - 2Btj + 2A = 0, (2.1) 

where A and B are constants of integration. The CH one-phase solution u{x, t) = 
(c — 2v) — 2ri{C,) is then obtained by inverting the differential of the third kind 

Jvo vw+WeTW^^W+^BT^) J no ^(e+i.)nLi(?-A:)' 

(2.2) 

where P1+P2+P3 = c—2u. Let —v < (3^ < < Pi, so that 77(C) is real periodic in the 

interval [03,(32]. Since — ^ has constant sign for 77 G [/33,/?2]i by a standard argument, 
a?7 

the (real) inverse function 77(C) exists and is monotone in C G [0, Z], where Z is the 
half period of the travelling wave solution. To invert (|2.2p and obtain 77 = 77(C) let us 
introduce the elliptic curve 

3 

£: {w' = B.(0 = {^ + 1^) - f^^)}' (2-3) 

1=1 

with homology basis of cycles a, b as in the figure below and let 



^OlW Ja W ^(/32 +J^)(/3l -/33) 

be the normalized holomorphic differential on f , with K{s) the complete Jacobi elliptic 
integral of modulus s 

2 ^ (/3i + ^)(/32-/33) , . 

ii32 + <^m-f33y 



Fig. 2.1. The homology basis 



Then inversion of the normahzed holomorphic differential 



is given by the Jacobi theorem [T^ and takes the form 



where 



lo dz 'i?3(z + po + i + 5) ■ 



K{s) 



(2.6) 



A(s,p)= / , P=^ , (2.7) 

io (l-pz2)^(l-z2)(l-s2z2) ^2+1^ 

i?3(z) is the third Jacobi-theta function defined by the Fourier series 

+00 K'{s] 
dziz.r) — 7 exp [jTrn^T + 27rmz] , t — i (2.8) 

^(s) 



if'(s) = K{'\/\ — s^) and po is defined by 

/■°° . 1 /■ . 1 /"V dz 

i-i. 2 Jp^ 2K(s) Jo ^(f - z2)(i _ 32^2) 

The function 77(2;) is an elhptic function in z with periods 1 and r and it is real periodic 
in z along the complex line z = t/2 + z' , z' E R. 
Finally, inserting (12.61) into (12. 2p . we get 



C = fc/o + ^)rfz ^ fc/o(t. + P3)^^^z - Hog ^! • (2-10) 
Jr/2 K{s) i?3(z+p,r) 

Equations (|2.6p and (|2.f Op give the one-phase solution u(C) in implicit form. The real 
half period, Z, of u{() is 

/•(r+l)/2 

Z = /c/o / (r/(z) + t/)(iz. (2.11) 

Jt/2 
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Normalizing the half period Z to tt, from (|2.1ip we obtain the wave number and the 
frequency 



27r 



{i + v)d\ \ _ ^v/(/32 + i^)(/3i-/33) 



/a /R(0 ; 2(/?3 + z^)A(s,p) ' (2.12) 

= (2i. + A + /32 + /?3)fc. 

Summarizing, the one-phase real nonsingular solution to the CH equation takes the 
form 

^(x,t,e)=/3i+fe+/33 + 2^-2(/?3 + ^)^^ + |-;^log ^^["~^'^j , 

e k d^{z+p,T) 

(2.13) 

Observe from the above formula that u{x,t,e) is not a meromorphic function of x 
and t. An alternative method to obtain u{x, t, e) is discussed in [3] where the solution 
is expressed in terms of conveniently generalized theta-functions in two variables, 
which are constrained to the generalized theta-divisor of a 2~dimensional generalized 
Jacobian. 

The envelope of the oscillations is obtained by the maximum and minimum values 
of the theta function and gives 

.™. = /3r+/3.+/33 + 2.-2(/33 + .^^^'^'^) ■ 



= /3i+/32 +/33 + 2z^-2(/?3 + z.) 



K{s) /o,?3(p+i)' 
A(,s,p) Ad'.jp) 
K{s) IoMp)' 



3. Small amplitude limit of the one-phase solution. We analyze the one- 
phase solution (|2.13p near the leading edge x^{t), namely when the oscillations go to 
zero. For this purpose we need to study the CH modulation equations more in detail. 

3.1. Camassa-Holm modulation equations. In [T], Abenda and Grava con- 
structed the one-phase CH modulation equations and showed that Pi, i = 1,2,3 are 
the Riemann invariants. The CH modulation equations take the Riemann invariant 
form 

a* A + C,(/3)9,A = 0, z = 1, . . . , 3, (3.1) 
where the speeds Ci(/3) in (|1.6p are explicitly given by the formula [T] 



Ci(/3) =/3i-f /32+A3 + 2i^ + 

C2{l3)=Pl+p2+P3+2l^ + 



{(32 + v)E{s) 
2(/32-/33)A(s,p) 



(/33 + y){(ii - 132) 



0*3 /3) = /?! 4- /?2 + /33 + 2i^ + 2—— — -— — ^rrr^, 

{(32 + i^)[K{s)- E{s)] 

with f3 = {f3i, (32, Pa), K{s), A(s, p) as before and E{s) the complete elliptic integral of 
the second kind with modulus s^. In the limit when two Riemann invariants coalesce, 
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the modulation equations reduce to the Hopf equation 

dtu + (3u + 2v)dxU = 0. 

The CH modulation equations are integrable via the generalized hodograph transform 
introduced by Tsarev [30j , 



x = ~C^{|3)t + w^{|3) i = l,2,3, 



(3.3) 



which gives the solution of (|3.ip in implicit form. The formula of the 'Wi{f3) for the 
specific CH case, has been obtained in [I], [T^ : 



wm = qif3) + (cm - + /32 + /33 + 2^.)) dpMl3) 
where the function q = q{/3) is given by 
1 



1,2,3, 



(3.4) 



(z(/3i,/32,/33) 



2V2J-1J-1 



with /_ the inverse of the decreasing part of the initial datum uo{x). The above 
formula of q is valid for x < XM{t) where XM{t) is the coordinate of the minimum 
value of the solution u{x,t) of the Hopf equation. For x > xuit) the corresponding 
formula for q contains also the increasing part of the initial datum (see [28j). In [19] 
it is shown that if /"'(uc) < 0, then the solution /3i(t) > /?2(i) > P^it) of the system 
(13.41) exists for some time t > tr- 



In order to take the small amplitude limit of the solution ()1.4p . we rewrite the 
system (|3.3|) in the equivalent form [19] 



{Cit + wi — x)a = 0, 

C2t + W2 — X = 0, 
{C2 - C3)t + {W2 - W3) 

P2-I33 



(3.6) 



where 



{P2 + y)E{s) 



Pi 



0, 



{f3i-(32){P3 + v)K{p,s) 
with s and p as in (jl.Sp and (j2.7p . and we perform the limit (5^0, where 

(i2=V + 6, ^ V - S, /?! = M. 

Let 



(3.7) 



Q{u,v) = q{u,v,v) 



V2 dp 

-1 Vl - M 



f- 



\ + p 1- p 
V H u 



$(m, v) = duQiu, v) + dyQ{u, w) = — 



V2 dp 



1 \/l - M 



l+p l-p ^ 
V H u 



(3.8) 
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then the foUowing identities hold 

f-{u) = Q{u,v) + 2{u - v)duQ{u,v), 
2duQ{u, v) - dyQ{u, v) 



duvQ{u,v) = 
duvvQ{u,v] 



2{u-v) ' (3.9) 

4:duvQiu,v) - dvvQiu,v) 



2{u-v) 

Substituting (|3.8p . (|3.9p and the expansion of the elhptic integrals E{s), K{s) and 
A(s, p) as s, p ^ [19 into (|3.6p . we arrive to the system 

. ^3. 4- 2.)t + /„(.) - ^ (3t + ^K^))(2^ + . + 3.) ^ ^ 

2 (u — v)^ 

U — V 

X =(3u + 2iy)t + /_(u) - 2 {t{u + 2v + + <^>{u, v){v + ly) + duQ{u,v){u - v)) 



{2v-u + v){it + $(u, v)) - 2{u -v){v + iy)dy^{u, v) 



2{u + v) 

{u ~v){v + v) - 8uv + Sviy - 2uv + iv"^ + "iu^ 



r (9™$(m,u)- 



2{u + v) 8{u + v){v + v){u ~ v) 



(3t + $(u,w)) +0{5-'), 



{2v-u + v){it + <^{u,v)) - 2{u - v){v + iy)dy^{u,v) ^ 



2{u + iy) 

(3.10) 

From the above, we deduce that, in the limit (5 — > 0, the hodograph transform (13. 3p 
reduces to the form 

= {3u{t) + 2v)t + /_ (u) - X, 

Q={u-v){t + duQiu, v)) + {v + z/)(3t + 3>(it, v)), (3.11) 
0={u-2v- i/)($(u, v) + 3i) + 2{v + t/)(u - v)dy^{u, v). 

The above system of equations enables one to determine x, u and v as functions of 
time. We denote this time dependence {t), V = v{t) and u = u{t). Observe 

that u{t) = u(x~ (t),t) where u{x,t) is the solution of the Hopf equation. In what 
follows we will always denote by u or u(t) the solution of the Hopf equation u{x, t) at 
the leading edge x = x~ (t) while we will refer to the solution of the CH equation as 
u{x, t, e). The derivatives with respect to time of these quantities are given by 

Xf =3u + 2v+ (3t + f'{u))ut, 
_ {u - v){u + 2v + 2,v) 

" ~ {u + v){it + f'_{u)) ' (3.12) 

A{u — — Av){v + z/)(u — v) 
Vt = y , 

with 

V = A{v + vf{u - vfdyy<^>{u, v) - {3t + $(w, v)){8viy + 3u^ ~ 8uv - 2uv + + Su^). 

(3.13) 

We are interested in studying the behaviour of the one-phase solution (|2.13p near 
the leading edge, namely when < a; — a;^ (i) <C 1. To this aim, we introduce two 
unknown functions of x and 

5 = 5{x^x~{t)), A = A(x- a;-(i)). 
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which tend to zero as x x (t). We now derive the dependence of A as a function 
of a; — x^{t). Let us fix 

(32 = V + 5, (33^v-6, /3i=u + A, A 0. (3.14) 

For the first equation in (|3.6p near 6 = 0, we find 

We then substitute fii = u + A in the above equation, we insert p.Sp . (|3.9p and (|3.1ip 
into it, and we get 



so that 



Similarly, using the second equation in (|3.6p . we arrive at 



A-^-^. (3.16) 

3i + /-(M) 



x-x-(t)«<52c, c=-— — — ^— ->0, (3.17) 

8(u + v)[v + z/)(m — v) 

with y defined in ((XT^ . 

Theorem 3.1. In the limit 

f32=v + S, f33 = v-6, /3i = u + A, (5, A 0, 

f/ie one-phase solution of the CH equation, implicitly defined by /las the following 

trigonometric expansion 

u(x,t,e) « u(t + -44 + 2(5 cos - - — 1 - cos 2- 

^ ' ^' 3t + fLiu) \ej 2 {u - v){v + ly) \ \e 

with 

e = eo + 6 + o((x~x_(t))2), 
Co - -2^{u-v)(v + jy)mu,v)+3t) a = 

\ V + ly 

Proof: We first prove ([XTO]) . Let = ^/e in pn]) . so that 

e = fc [x - (/3i + /?2 + /33 + 2;/)< - g(/3i, /32, /Sg)] , 

with k as is (12.41) . g as in (13. 5p . In the limit S,A^O, with f]2 = v + S, (3s=v~5, 
/3i = u + A, we get 



(3.18) 



(3.19) 



ti + 1/ \ 2(u — t^) / \ 16 (w + i^)^(M — v)"^ 
'{u + 2iy + 2v + A)t - Q{u, v) - A(9„Q(u, v) - — 9™Q(u, w) ) + 0{S'^), 
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where we use A — 0{d'^) to estimate the error in the formula above. Then, we insert 
the solution x_{t) of equations (|3.1ip . we use p.l6p and (|3.17p to estimate A and 6^ 
and we get 

e = -2^{u-v)iv + iy){<^{u, v) + 3t) + - x_ {t)) + 0((x - x_ {t)f), 

y V + 1' 

where $ is defined in (|3.8p . which coincides with p.lPp . 

We now prove that in the limit 6 —^ 0, u{x,t,e), defined implicitly as a function 
of C in p.4p . has the following trigonometric series expansion 



u 



{x, t, 6) = /3i + 2<5cos(C) - ^ . (1 - cos(2C)) + 0{5'). (3.20) 

2 - v){v + v) 

To this aim, we expand u{x,t,e) and C defined in (|1.4p in trigonometric series of z 
and we compute z as a function of 

We need the following series expansions for the Jacobi theta functions [32] 



(3,21) 



^ ' n— 1 

where j> is as in (|2.9p and 

C = exp (-.i.'(,s)/i.(.)) ^ 3(^1%;^^) ^ + (3.22) 
Inserting the above expansion and 



sm(27rp) = 2 — h 0((5 ), 

Pi + 



sin(4.,) = -4 ^^- + -^;y(\t^ ^^^ + 0(.^), 
into (|3.2ip . we get 

u{x,t, e) = /3i + 25 cos(2^z) + ^ J?L_il_^ (1 - cos(4^z)) + 0{S'), 
C = 27rz sin(27rz) + 0(6^) 



from which we find 



2^^ = C+"-^S + 0{S% 

V + v 



and the assertion easily follows. □ 
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4. Painleve equations at the leading edge. In this section we propose a 
multiscale description of the osciUatory behavior of the solution to the CH equation 
in the small dispersion limit (e — > 0) close to the leading edge X- {t) where (32 = 03 = v 
and /?i = u. We follow closely the approach [18] for the corresponding KdV situation. 
The ansatz for the multiscale expansion to the CH solution close to the leading edge 
is inspired by the asymptotic solution in the Whitham zone discussed in the previous 
section. Numerically we find that the quantity 5 in (|3.18p is of the order e-^^^. This 
implies with (|3.17p that x — x^{t) « e^^^ in the double scaling limit e ^ and 
X — » x-{t). We are thus led to introduce the rescaled coordinate y near the leading 
edge, 

y^e-'/^{x-x^(t)), (4.1) 
which transforms the CH equation (|l.ip to the form 

(3U + 2l/ — X^)Uy + e^^^{Ut — 2UyUyy — (m — X^)Uyyy) " E^^^Uyyt — 0, (4.2) 

where x^ = ■^x^{t). 

Numerically the corrections to the Hopf solution near the leading edge are of 
order e^/'^ and thus we make as in the ansatz 

u{y, t, e) = C/o(y, t) + e^'^U^{y, t) + e^"'U2{y, t) + eU^iy, t) + . . . , (4.3) 

where Uo = u{t) is the solution of the Hopf equation at the leading edge. We assume 
that the terms Uk, k > 1 contain oscillatory terms of the order 1/e and take the form 



U,iy,t)^aiy,t)cos(^^i^ 
U2iy, t) = hi{y, t) + 62(2/, t) cos 



2V'(2/,i) 



e 



Uziy, t) = co(j/, t) + ci (y, t) cos [t^^] + 02(2/, t) sin ' 



(4.4) 



+ cz{y,t) cos 

where 



iPiy, t) = ^0 (y, t) + e^/Vi (y, t) + e'/V2(y, t) + £^'3(2/, t) + ■ ■ ■ ■ (4.5) 

Terms proportional to sin(?/'/e) in \J\ can be absorbed by a redefinition of ?/;, the 
higher order terms are chosen to compensate terms of lower order appearing in the 
solution to CH due to the non-linearities. We immediately find c^y'i/'o = = dyili\ 
from the terms of order and . The term of order ^ gives 

(1 + (5y^2)') dt^^ + dyi)2 (3u + 2v~ x-^) + {u{t) - X-,) {dy^2f = 0. (4.6) 

In particular, (j4.6p is algebraic in 9j^?/'2 with coefficients depending on t only, i.e., 
ip2{y,t) = f2it)y + /i(i), with fi.2{t) to be determined. At the order e^/'^, we get 
dtipi =0 and 



a'(l + /2)'-8fe2/|(« + ^), 



(4.7) 
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f^{x^ -u) + 2/2 (xr -2v- -iu = 0. (4.8) 

At the order e^^'^, using (|4.6p to (14. 8p . we find dy^z — and 



adya (/I - 1) _ a3 (3 ^ 7^2) + ^. 



,3 



l («2(3 + /|)+4u,,y)(l + /|)2 
' "8 /|(3 + /!)(. + .) + (^-^^^ 



with integration constant, 



„ _ + 3 , ay (1 + /I) 



(4.11) 



= 0. 



32(u + ^.)2/r " 4 /|(-3 + /|)(u + z.)2) 
g (1 + /|)39,V,2 ffia (3 + /|)(l + /|)2 
2/2(u + i/)(-3 + /|) 2 (7. + j/)(-3 + /|) 

We determine 9j,'(/'2 inserting (|3.12p inside (|4.8|) . so that 

^yV'l = /I e ^ , 4.12 

[v + iyu— u — 2v} 

where the sign in front of the square root has to be chosen in such a way that the 
r.h.s. is positive. Then comparing (|4.12p . (14. 5p . p.l9p and (|4.ip . we conclude that 



My,t)^yJ^. (4.13) 



Inserting ff^J^ into (gjl]), we find 



5.^0(0 = -4./^^. (4.14) 

V u — 71 + ;/ 

By comparison with (j3.19p . we observe that 

Then integrating the l.h.s. of the above expression between tc and t we arrive at (j3.19p . 
namely 



Mt) = - v)iv + iy){^{u,v) + 3t). (4.15) 

Moreover, consistency between (|4.3p and (|3.18p implies that 



5=1^, (4.16) 



from which we immediately verify that 

2/3, u + 



2 (w-u)(t; + j/)' 
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3t + f'{u) 2 {u~v){v + i^)' 



with .91 (t) = in (|4T0)) . 
Summarizing, we get 

u{x,t,e) = u(<)+e^/^acos (-]+e'^^^ 



a'^{u + v){cos(2'ip/e) - 1) 



+0(e), 
(4.17) 



where 



^{y, t) = -2^iu-v)iv + iy)mu, v) + it) + e^/^y^"^ + 0(e), (4.18) 



with Phi{u,v) defined in (|3.8p and a satisfies the Painleve-II equation 



dyya — Aiay — A2a'^ — 0, 



with 



{u + vfvt 



{u + uf 



A(u - v){v + i/)3(3j/ + 4i; - u) 
Using (|3.12p . we get 

{u + vY 



'i2{u - v){v + vf 



(4.19) 



(4.20) 



Ai - 



y(w + !/)2' 

with V defined in p.l3p . Making the substitution z = ay, a = j3A with 

a^/2 



a={A^f\ (3 



we arrive to the special Painleve II equation in normal form 



(4.21) 



For X — x^ify > 0, from the small amplitude limit of the one-phase solution to the 
CH equation we get 



with c defined in (I3.17p . so that in the limit e 0+ (equivalently y — > +oo or 
z — > — oo), we have 



Aiz) » ^-z/2. 

If X ^ X-{t), the CH solution is approximately the solution to the Hopf equation so 
that a{y) ~ 0, for y — > —oo (equivalently z — > +oo). Therefore we conclude that 

lim A{z) = 0. 

z — * + oo 

The solution to the Painleve II equation satisfying such asymptotic condition exists 
and is unique and pole free on the real line according to Hastings and McLeod [3T]. 
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5. Numerical solution of CH and Whitham equations. In this section 
we will solve numerically the CH and the Whitham equations for initial data in the 
Schwartzian class of rapidly decreasing functions with a single negative hump. As a 
concrete example we will study the initial datum 

uq :— u{x,0) = — sech^a;. (5-1) 

For this initial datum the non-breaking condition (|1.2p is satisfied if < (i^ — l)/2- 

5.1. Numerical solution of the CH equation. The resolution of the rapid 
modulated oscillations in the region of a dispersive shock is numerically demanding. 
The strong gradients in the oscillatory regions require efficient approximation schemes 
which do not introduce an artificial numerical dissipation into the system. We there- 
fore use Fourier spectral methods which are known for their excellent approximation 
properties for smooth functions whilst minimizing the introduction of numerical vis- 
cosity. We restrict to initial data, where uq — e^uo^xx + v does not change sign, to 
ensure analyticity of the CH solutions. The Schwartzian solutions can be treated as 
effectively periodic if the computational domain is taken large enough that the solu- 
tion is of the order of machine precision 10"^^ in Matlab) at the boundaries. We 
always choose the computational domain in this way to avoid Gibbs phenomena at 
the boundaries. 

Even with spectral methods, a large number of Fourier modes is needed to resolve 
the rapid oscillations numerically. To obtain also a high resolution in time, we use 
high order finite difference methods, here a fourth order Runge-Kutta scheme. For 
stability reasons a sufficiently small time step has to be chosen. Unconditionally stable 
implicit schemes could be used instead, but these can be only of second order. Such 
approaches would be too inefficient for the precision requested. Notice that the terms 
with the highest derivative in CH p.l|l are not linear in contrast to the KdV equation 
(we choose the KdV equation in a way that it has the same dispersionless equation 
as (jl.ip . the term proportional to v can be always eliminated here in contrast to CH 
by a Galilean transformation) 

ut + (3u + 2v)ux -I- <?Uxxx = 0. (5.2) 

For such equations efficient integration schemes are known. In |24| it was shown that 
exponential time differentiation (ETD) methods [H I26j are the most efficient in the 
small dispersion limit of KdV. For CH such an approach is not possible due to the 
nonlinearity of the highest order derivatives. But we can use analytic knowledge of 
the solution: The dispersionless equation, the Hopf equation, will have a point of 
gradient catastrophe at the critical time tc, for the example (|5.1[) tc = •\/3/4 « 0.433. 
For times t <^ tc, the CH solution is very close to the Hopf solution with moderate 
gradients. In this case we can use a larger time step. For times close to the critical 
time and beyond (in our example t > 0.4), we have to use considerably smaller time 
steps. 

The main difference to the KdV equation is the non-locality of CH, which provides 
some filtering for the high frequencies: if we denote the Fourier transform of u with 
respect to x by u{k,t) := u{x,t)e^'^^dx, we get for (|l.ip in Fourier space 

= I ^ g2fc2 + 2mku - €^{2u^u^x + S^w) j • 

The term 1/(1 + e^A:^) in the above equation is the reason why CH gives for large 
spatial frequencies a better approximation to one-dimensional wave phenomena than 
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KdV with a term proportional to in Fourier space. Its effect in the present context 
is twofold: First it provides a high frequency filtering which allows in practice for 
larger time steps in the computation. Secondly it suppresses the rapid modulated 
oscillations in the shock region of the dispersionless equation. This can be seen in 
Fig. 15. li where solutions to KdV and CH for the same initial data are shown. The 

CH, negative data CH, positive data 




-5 5 -4 -2 2 4 6 



KdV, negative data KdV, positive data 




-5 5 -4 -2 2 4 6 



X X 

Fig. 5.1. Solution to the CH equation il.lV and KdV equation i5. 2V with v = 1.2 for the initial 
data uq = —sech^x on the left side and uq = sech^x on the right side at t = 1 for e = 10~^. 

KdV solutions always show much more oscillations than the CH solutions. 

Notice that despite the lower number of oscillations in the CH solutions in Fig. 15.11 
depending on the value oi v & higher number of Fourier modes is needed to resolve 
the oscillations. This is due to the fact that the CH solution is less smooth than 
the KdV solution and thus less localized in Fourier space, as can be seen in Fig. 15.21 
Therefore, to treat the same values of the small dispersion parameter e for CH as 
for KdV, we need higher temporal and spatial resolution. It is thus computationally 
more demanding to obtain the same number of oscillations in CH solutions as in KdV 
solutions, the latter being important to obtain a valid statistics for the scaling studied 
below. 

The quality of the numerics is controlled via energy conservation for CH, 

E 2vu^ + e^uul + u^. 

The numerically computed energy will depend on time due to numerical errors. As 
was discussed in [24) . energy conservation can thus be used to check numerical accu- 
racy. In practice energy conservation overestimates numerical precision by 1-2 order 
of magnitude. We typically solve the CH equation with a relative numerical error 
/S.E / E = 10^^. This ensures that the difference between the numerical and the 
asymptotic CH solution, which is typically of the order of e or larger, is entirely due 
to the asymptotic description. 

5.2. Numerical solution to the Whitham equations. The Whitham equa- 
tions p.6p have a similar form as the respective equations for KdV. We use the same 
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CH, negative data 




-0.8 -0.6 -0.4 -0.2 



KdV, negative data 



CH, positive data 




KdV, positive data 



3.5 4 4.5 



Fig. 5.2. Magnification of the oscillations in Fig. 15. l\ 



procedure to solve them numerically: we first solve the equations at the edges of the 
Whitham zone and then at intermediate points. For details the reader is referred to 
[T8] . Typical solutions for the Whitham equations can be seen in Fig. 15.31 In the 
shown example, the quantity f]^ crosses the hump at —1. 




-0.6 -0.5 



-0.2 -0.1 




Fig. 5.3. Solutions /3i, i = 1,2,3 to the Whitham equations Sl.dt for 1x0(3;) = —sech?x and 
t = 1 for two values of v. 

In contrast to KdV, the system (|1.6p is not strictly hyperbolic, i.e., the speeds C^, 
i — 1,2,3, in (|3.2p do not satisfy for all x and t > tc the relation Ci > C2 > C3. In 
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fact the lines Ci {x) can cross for a given time for v ^ 1 as can be seen in Fig. 15.41 For 




Fig. 5.4. The speeds d, « = 1,2,3 13. S\) for uq{x) = —sech x and t = 1 for two values of v. 

V = 1.02 the speeds Ci and C3 intersect in the shown example in the interior of the 
Whitham zone as can be seen in more detail in Fig. 15.51 This behavior of the speeds 




Fig. 5.5. Magnification of the lines in the lower figure of Fig. \ 5.4\ close to the crossing of the 
lines, 

has no influence on the numerical solubility of the Whitham equations. It also does 
not influence the quality of the asymptotic solution in these cases as can be seen in 
Fig. EH 

5.3. Quantitative comparison of the CH solution and the asymptotic 
solution. The asymptotic description of the small dispersion limit of the CH equation 
is as follows: for times t < tc, the Hopf solution for the same initial datum provides 
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an asymptotic description. For t > t^, the Whitham zone opens. Outside this zone, 
the Hopf solution again serves as an asymptotic solution. In the interior the one- 
phase solution to the CH equation describes the oscillatory behavior. It is given on 
an elliptic curve with branch points being solutions of the Whitham equations. 

Below we will study the validity of this asymptotic description in various regions 
of the (x, t)-plane. To study the e dependence of a certain quantity A, we perform a 
linear regression analysis for the dependence of the logarithms, Inyl = alne + 6. We 
compute all studied quantities for the e values e = 10^" with a G [1, 1.25, 1.5, . . . , 3]. 
Generally it is found that the correlations and the standard deviations are worse than 
in the KdV case due to the lower number of oscillations. 

Before breakup, t <^ tc : 
For times much smaller than the critical time, we find that the Loo norm of the 
difference between Hopf and CH solution decreases as e^. More precisely we find 
by linear regression an exponent a = 1.91 with correlation coefficient r — 0.999 and 
standard deviation aa = 0.06. 

At breakup, t ^ tc : 

For times close to the breakup time, the Hopf solution develops a gradient catastrophe. 
The largest difference between Hopf and CH solution can be found close to the breakup 
point. We determine the scaling of the L^o norm of the difference between Hopf 
and CH solution on the whole interval of computation. We find that its scaling 
is compatible with e^^^ as conjectured in [TT]. More precisely we find in a linear 
regression analysis a = 0.28 (2/7 — 0.2857 . . .) with a correlation coefficient r = 0.998 
and standard deviation aa — 0.015. An enhanced asymptotic description of the CH 
solution near the breakup point in terms of a solution to the Painleve 12 equation was 
conjectured in [HI and studied numerically in |17| . 
Times t > tc'- 

For times t ^ tc it can be seen from Fig. ll.2l and Fig. l5.6l that the asymptotic solution 
gives a very satisfactory description of the oscillations except at the boundaries of the 
Whitham zone. The asymptotic solution is so close to the approached solution that 
one can only see discrepancies near the boundary of the Whitham zone, where the 
asymptotic solution is just C°. Thus one has to consider the difference between the 
solutions as shown in Fig. 15.71 The quality of the numerics allows the study of the 
scaling behavior at various points in the Whitham zone. 

Interior of the Whitham zone: 
If we study the e-dependence of the Loo-norm of the difference near the middle of the 
Whitham zone (we take the maximum of this difference near the geometric midpoint), 
we find that the norm scales as e. More precisely we find an exponent a — 0.9 with 
correlation coefficient r = 0.999 and standard deviation aa — 0.04. There is obviously 
a certain arbitrariness in our definition of this error especially for larger values of e, 
where there are only few oscillations. In these cases the errors are read off close to 
the boundaries of the Whitham zone where much smaller exponents for the error are 
observed (see below). This scaling gives nonetheless strong support for the conjectured 
form of the phase of the one-phase solution in the Whitham zone, since even a small 
analytical error in the phase would lead to large errors which would not decrease with 
e. 

Leading edge of the Whitham zone: 
Oscillations can always be found outside the Whitham zone, whereas the Hopf solution 
does not show any oscillations. The biggest difference always occurs at the boundary 
of the Whitham zone. It scales as e^/^. More precisely we find in the Hopf zone 
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v=1 .2 




-1 I 1 1 1 1 1 1 ^ 

-0.9 -0.8 -0.7 -0.6 -0.5 -0.4 -0.3 

X 



Fig. 5.6. CH solution (blue) and asymptotic solution (magenta) for the initial datum uq = 
— sech^x and e = 10~^ for t = 1 and two values of v 




-0.1 -0.1 



-0.2' • • ' -0.2' • ^ 

-1 -0.5 -1 -0.5 

X X 



Fig. 5.7. Difference of the CH solution and the asymptotic solution for the initial datum 
UQ = —sech^x and u = 1.2 for t = 1 for several values of e. In the Whitham zone (red), the 
asymptotic solution is given in terms of the one-phase CH solution solution, outside by the solution 
to the Hopf equation with the same initial datum. 



a = 0.33 with r = 0.99 and aa = 0.04. In the interior of the Whitham zone, a has 
a similar value, but the correlation is worse. If one studies the scaling of the zone, 
where the difference between Hopf and CH solution near the leading edge has absolute 
value larger than some treshold (we use 10"**), we find a decrease compatible with 
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g2/3^ more precisely a = 0.81 with r — 0.99 and = 0.07. 

Trailing edge of the Whitham zone: 
The biggest difference is always found at the boundary of the Whitham zone. Its 
scaling is compatible with e^/^. We find a — 0.51 with r — 0.9997 and (Tq — 0.01. 

6. Numerical study of the multiscales expansion for the CH equation. 

In this section we study numerically the multiscales solution to CH derived in sec- 
tion [D It will be shown that the latter provides a better description of the asymptotic 
behavior near the leading edge of the Whitham zone than the Hopf or the one-phase 
solution to CH as discussed in section [5l We identify the zone, where the latter gives 
a better description of CH than the former and study the e-dependence of the errors. 
For the numerical computation of the Hastings-McLeod solution and a comparison to 
the KdV case, we refer the reader to [T5j . 

In Fig. 16.11 we show the CH solution, the asymptotic solution via Whitham and 
Hopf and the multiscales solution near the leading edge of the Whitham zone. It can 
be seen that the one-phase solution gives a very good description in the interior of 
the Whitham zone as discussed in section [5l whereas the multiscales solution gives as 
expected a better description near the leading edge. 



1 




Fig. 6.1. The figure shows in the upper part the numerical solution to the CH equation for 
the initial datum uq = —sech}x and v = 1.2, e = 10~^ at t = 1, in the middle the corresponding 
asymptotic solution in terms of Hopf and one-phase solution, and in the lower part the multiscales 
solution. 

In Fig. l6.2l the CH and the multiscales solution are shown in one plot for e — 10^^. 
It can be seen that the agreement near the edge of the Whitham zone is so good that 
one has to study the difference of the solutions. The solution only gives locally an 
asymptotic description and is quickly out of phase for larger distances from the leading 
edge. The difference between CH and multiscales solution is shown for several values 
of e in Fig. 16.31 It can be seen that the maximal error still occurs close to the Whitham 
edge, but that it decreases much faster with e than the error given by the Hopf and 
the Whitham solution. A linear regression analysis for the logarithm of the difference 
A between CH and multiscales solution near the edge gives a scaling of the form 
A oc e° with a — 0.58 with standard deviation aa = 0.17. Since there are much less 
oscillations in the CH case than in KdV, the found statistics is considerably worse 
in the former case than in the latter |18| . which is reflected by the low correlation 
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Fig. 6.2. The numerical solution to the CH equation for the initial datum uq = —sech}x and 
V = 1.2, e = 10"'^ at t = 1 in blue and the corresponding multiscales solution in green. 

coefficient r ~ 0.93 and the comparatively large standard deviation. Nonetheless the 
found scaling is in accordance with the e^/'^ scaling expected from the multiscales 
expansion. 




-0.8 -0.7 -0.6 -0.8 -0.7 -0.6 
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-0.8 -0.7 -0.6 -0.8 -0.7 -0.6 



Fig. 6.3. The difference between the numerical solution to the CH equation for the initial datum 
uq = —sech^x and v = 1.2 at t = 1 and the corresponding multiscales solution for 4 values of e. The 
interior of the Whitham zone is shown in green, the exterior in blue. 

As can be already seen from Fig. 16.11 the multiscales solution gives a better 
asymptotic description of CH near the leading edge of the Whitham zone than the 
Hopf and the one-phase solution. This is even more obvious in Fig. 16.41 where the 
difference between CH and the asymptotic solutions is shown. This suggests to identify 
the regions where each of the asymptotic solutions gives a better description of CH 
than the other. The results of this analysis can be seen in Fig. 16.51 This matching 
procedure clearly improves the CH description near the leading edge. In Fig. 16.61 we 
see the difference between this matched asymptotic solution and the CH solution for 
two values of e. Visibly the zone, where the solutions are matched, decreases with e. 
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Fig. 6.4. The difference between the numerical solution to the CH equation for the initial 
datum uo = —sech^x and v = 1.2, e = 10~^ at t = 1 for t = 10~^ and the corresponding multiscales 
solution in blue, and the difference between CH and Hopf and one-phase solution in green. 




-1 -0.8 -0.6 -0.4 -0.2 0.2 0.4 




0.2 0.4 



Fig. 6.5. In the upper part one can see the difference between the numerical solution to the CH 
tion for the initial datum uq = —sech^x and v = 1.2, e = lO"'^ at t = 1 and the corresponding 
asymptotic solution in terms of Hopf and one-phase solution. The lower figure shows the same 
difference, which is replaced close to the leading edge of the Whitham zone by the difference between 
CH solution and multiscales solution (shown in red where the error is smaller than the one shown 
above). 



There is a certain ambiguity in the precise definition of this matching zone due 
to the oscihatory character of the solutions. Because of the much lower number of 
oscillations than in KdV, the statistics is considerably worse in the CH case than in 
the KdV case. The limits of the matching zone for several values of e can be seen in 
Fig. 16.71 Due to the lower number of oscillations in the Hopf region, the matching 
zone extends much further into this region than in the Whitham region. The width of 
this zone scales like with a = 0.51 and standard deviation <Ja — 0.06 and correlation 
coefficient r = 0.99. 
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Fig. 6.6. The difference between the numerical solution to the CH equation j or the initial datum 
uq = —sech^x and v = 1.2 at t = 1 and the corresponding asymptotic solution in terms of Hop f and 
one-phase solution in blue and CH and multiscales solution in red, where the latter error is smaller 
than the former, for two values of e . 
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Fig. 6.7. The edges of the zone where the multiscales solutions provides a better asymptotic 
description of CH than the Hopf or the one-phase solution for the initial datum uq = —sech^x and 
V = 1.2 at t = I in dependence of e. 
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